home *** CD-ROM | disk | FTP | other *** search
/ EnigmA Amiga Run 1999 March / EnigmA AMIGA RUN 35 (1999)(G.R. Edizioni)(IT)[!][issue 1999-03].iso / earcd / grafica / pvrgjpeg / transform.c < prev   
C/C++ Source or Header  |  1999-01-01  |  12KB  |  525 lines
  1. /*************************************************************
  2. Copyright (C) 1990, 1991, 1993 Andy C. Hung, all rights reserved.
  3. PUBLIC DOMAIN LICENSE: Stanford University Portable Video Research
  4. Group. If you use this software, you agree to the following: This
  5. program package is purely experimental, and is licensed "as is".
  6. Permission is granted to use, modify, and distribute this program
  7. without charge for any purpose, provided this license/ disclaimer
  8. notice appears in the copies.  No warranty or maintenance is given,
  9. either expressed or implied.  In no event shall the author(s) be
  10. liable to you or a third party for any special, incidental,
  11. consequential, or other damages, arising out of the use or inability
  12. to use the program for any purpose (or the loss of data), even if we
  13. have been advised of such possibilities.  Any public reference or
  14. advertisement of this source code should refer to it as the Portable
  15. Video Research Group (PVRG) code, and not by any author(s) (or
  16. Stanford University) name.
  17. *************************************************************/
  18.  
  19. /*
  20. ************************************************************
  21. transform.c
  22.  
  23. This file contains the reference DCT, the zig-zag and quantization
  24. algorithms.
  25.  
  26. ************************************************************
  27. */
  28.  
  29. /*LABEL transform.c */
  30.  
  31. /* Include files */
  32.  
  33. #include "globals.h"
  34. #include "dct.h"
  35. #include <math.h>
  36.  
  37. /*PUBLIC*/
  38.  
  39. extern void ReferenceDct();
  40. extern void ReferenceIDct();
  41. extern void TransposeMatrix();
  42. extern void Quantize();
  43. extern void IQuantize();
  44. extern void PreshiftDctMatrix();
  45. extern void PostshiftIDctMatrix();
  46. extern void BoundDctMatrix();
  47. extern void BoundIDctMatrix();
  48. extern void ZigzagMatrix();
  49. extern void IZigzagMatrix();
  50. extern int *ScaleMatrix();
  51. extern void PrintMatrix();
  52. extern void ClearMatrix();
  53.  
  54. static void DoubleReferenceDct1D();
  55. static void DoubleReferenceIDct1D();
  56. static void DoubleTransposeMatrix();
  57.  
  58. /*PRIVATE*/
  59.  
  60. /* The transposition indices */
  61.  
  62. int transpose_index[] =          /* Is a transpose map for matrix transp. */
  63. {0,  8, 16, 24, 32, 40, 48, 56,
  64.  1,  9, 17, 25, 33, 41, 49, 57,
  65.  2, 10, 18, 26, 34, 42, 50, 58,
  66.  3, 11, 19, 27, 35, 43, 51, 59,
  67.  4, 12, 20, 28, 36, 44, 52, 60, 
  68.  5, 13, 21, 29, 37, 45, 53, 61,
  69.  6, 14, 22, 30, 38, 46, 54, 62,
  70.  7, 15, 23, 31, 39, 47, 55, 63};
  71.  
  72. int zigzag_index[] =              /* Is zig-zag map for matrix -> scan array */
  73. {0,  1,  5,  6, 14, 15, 27, 28,
  74.  2,  4,  7, 13, 16, 26, 29, 42,
  75.  3,  8, 12, 17, 25, 30, 41, 43,
  76.  9, 11, 18, 24, 31, 40, 44, 53,
  77. 10, 19, 23, 32, 39, 45, 52, 54,
  78. 20, 22, 33, 38, 46, 51, 55, 60,
  79. 21, 34, 37, 47, 50, 56, 59, 61,
  80. 35, 36, 48, 49, 57, 58, 62, 63};
  81.  
  82. int izigzag_index[] =
  83. {0,  1,  8, 16,  9,  2,  3, 10,
  84. 17, 24, 32, 25, 18, 11,  4,  5,
  85. 12, 19, 26, 33, 40, 48, 41, 34,
  86. 27, 20, 13,  6,  7, 14, 21, 28,
  87. 35, 42, 49, 56, 57, 50, 43, 36,
  88. 29, 22, 15, 23, 30, 37, 44, 51,
  89. 58, 59, 52, 45, 38, 31, 39, 46,
  90. 53, 60, 61, 54, 47, 55, 62, 63};
  91.  
  92. /*Some definitions */
  93.  
  94. #define MakeMatrix() (int *) calloc(BLOCKSIZE,sizeof(int))
  95. #define FixedMultiply(s,x,y)  x = ((x * y) >> s);
  96. #define DCT_OFFSET 128
  97.  
  98. /*START*/
  99.  
  100. /*BFUNC
  101.  
  102. ReferenceDct() does a reference DCT on the input (matrix) and output
  103. (new matrix).
  104.  
  105. EFUNC*/
  106.  
  107. void ReferenceDct(matrix,newmatrix)
  108.      int *matrix;
  109.      int *newmatrix;
  110. {
  111.   BEGIN("ReferenceDct");
  112.   int *mptr;
  113.   double *sptr,*dptr;
  114.   double sourcematrix[BLOCKSIZE],destmatrix[BLOCKSIZE];
  115.  
  116.   for(sptr=sourcematrix,mptr=matrix;mptr<matrix+BLOCKSIZE;mptr++)
  117.     {                             /* Convert to doubles */
  118.       *(sptr++) = (double) *mptr;
  119.     }
  120.   for(dptr = destmatrix,sptr=sourcematrix;
  121.       sptr<sourcematrix+BLOCKSIZE;sptr+=BLOCKWIDTH)
  122.     {                             /* Do DCT on rows */
  123.       DoubleReferenceDct1D(sptr,dptr);
  124.       dptr+=BLOCKWIDTH;
  125.     }
  126.   DoubleTransposeMatrix(destmatrix,sourcematrix);  /* Transpose */
  127.   for(dptr = destmatrix,sptr=sourcematrix;
  128.       sptr<sourcematrix+BLOCKSIZE;sptr+=BLOCKWIDTH)
  129.     {                             /* Do DCT on columns */
  130.       DoubleReferenceDct1D(sptr,dptr);
  131.       dptr+=BLOCKWIDTH;
  132.     }
  133.   DoubleTransposeMatrix(destmatrix,sourcematrix);  /* Transpose */
  134.   for(sptr = sourcematrix,mptr=newmatrix;
  135.       mptr<newmatrix+BLOCKSIZE;sptr++)
  136.     {    /* NB: Inversion on counter */
  137.       *(mptr++) = (int) (*sptr > 0 ? (*(sptr)+0.5):(*(sptr)-0.5));
  138.     }
  139. }
  140.  
  141. /*BFUNC
  142.  
  143. DoubleReferenceDCT1D() does a 8 point dct on an array of double
  144. input and places the result in a double output.
  145.  
  146. EFUNC*/
  147.  
  148. static void DoubleReferenceDct1D(ivect,ovect)
  149.      double *ivect;
  150.      double *ovect;
  151. {
  152.   BEGIN("DoubleReferenceDct1D");
  153.   double *mptr,*iptr,*optr;
  154.  
  155.   for(mptr=DctMatrix,optr=ovect;optr<ovect+BLOCKWIDTH;optr++)
  156.     {                           /* 1d dct is just matrix multiply */
  157.       for(*optr=0,iptr=ivect;iptr<ivect+BLOCKWIDTH;iptr++)
  158.     {
  159.       *optr += *iptr*(*(mptr++));
  160.     }
  161.     }
  162. }
  163.  
  164. /*BFUNC
  165.  
  166. ReferenceIDct() is used to perform a reference 8x8 inverse dct.  It is
  167. a balanced IDCT. It takes the input (matrix) and puts it into the
  168. output (newmatrix).
  169.  
  170. EFUNC*/
  171.  
  172. void ReferenceIDct(matrix,newmatrix)
  173.      int *matrix;
  174.      int *newmatrix;
  175. {
  176.   BEGIN("ReferenceIDct");
  177.   int *mptr;
  178.   double *sptr,*dptr;
  179.   double sourcematrix[BLOCKSIZE],destmatrix[BLOCKSIZE];
  180.  
  181.   for(sptr = sourcematrix,mptr=matrix;mptr<matrix+BLOCKSIZE;mptr++)
  182.     {
  183.       *(sptr++) = (double) *mptr;
  184.     }
  185.   for(dptr = destmatrix,sptr=sourcematrix;
  186.       sptr<sourcematrix+BLOCKSIZE;sptr+=BLOCKWIDTH)
  187.     {
  188.       DoubleReferenceIDct1D(sptr,dptr);
  189.       dptr+=BLOCKWIDTH;
  190.     }
  191.   DoubleTransposeMatrix(destmatrix,sourcematrix);
  192.   for(dptr = destmatrix,sptr=sourcematrix;
  193.       sptr<sourcematrix+BLOCKSIZE;sptr+=BLOCKWIDTH)
  194.     {
  195.       DoubleReferenceIDct1D(sptr,dptr);
  196.       dptr+=BLOCKWIDTH;
  197.     }
  198.   DoubleTransposeMatrix(destmatrix,sourcematrix);
  199.   for(sptr = sourcematrix,mptr=newmatrix;mptr<newmatrix+BLOCKSIZE;sptr++)
  200.     {    /* NB: Inversion on counter */
  201.       *(mptr++) = (int) (*sptr > 0 ? (*(sptr)+0.5):(*(sptr)-0.5));
  202.     }
  203. }
  204.  
  205. /*BFUNC
  206.  
  207. DoubleReferenceIDct1D() does an 8 point inverse dct on ivect and
  208. puts the output in ovect.
  209.  
  210. EFUNC*/
  211.  
  212. static void DoubleReferenceIDct1D(ivect,ovect)
  213.      double *ivect;
  214.      double *ovect;
  215. {
  216.   BEGIN("DoubleReferenceIDct1D");
  217.   double *mptr,*iptr,*optr;
  218.  
  219.   for(mptr = IDctMatrix,optr=ovect;optr<ovect+BLOCKWIDTH;optr++)
  220.     {
  221.       for(*optr=0,iptr=ivect;iptr<ivect+BLOCKWIDTH;iptr++)
  222.     {
  223.       *optr += *iptr*(*(mptr++));
  224.     }
  225.     }
  226. }
  227.  
  228. /*BFUNC
  229.  
  230. TransposeMatrix transposes an input matrix and puts the output in
  231. newmatrix.
  232.  
  233. EFUNC*/
  234.  
  235. void TransposeMatrix(matrix,newmatrix)
  236.      int *matrix;
  237.      int *newmatrix;
  238. {
  239.   BEGIN("TransposeMatrix");
  240.   int *tptr;
  241.  
  242.   for(tptr=transpose_index;tptr<transpose_index+BLOCKSIZE;tptr++)
  243.     {
  244.       *(newmatrix++) = matrix[*tptr];
  245.     }
  246. }
  247.  
  248. /*BFUNC
  249.  
  250. DoubleTransposeMatrix transposes a double input matrix and puts the
  251. double output in newmatrix.
  252.  
  253. EFUNC*/
  254.  
  255. static void DoubleTransposeMatrix(matrix,newmatrix)
  256.      double *matrix;
  257.      double *newmatrix;
  258. {
  259.   BEGIN("DoubleTransposeMatrix");
  260.   int *tptr;
  261.  
  262.   for(tptr=transpose_index;tptr<transpose_index+BLOCKSIZE;tptr++)
  263.     {
  264.       *(newmatrix++) = matrix[*tptr];
  265.     }
  266. }  
  267.  
  268. /*BFUNC
  269.  
  270. Quantize() quantizes an input matrix and puts the output in qmatrix.
  271.  
  272. EFUNC*/
  273.  
  274. void Quantize(matrix,qmatrix)
  275.      int *matrix;
  276.      int *qmatrix;
  277. {
  278.   BEGIN("Quantize");
  279.   int *mptr;
  280.  
  281.   if (!qmatrix)
  282.     {
  283.       WHEREAMI();
  284.       printf("No quantization matrix specified!\n");
  285.       exit(ERROR_BOUNDS);
  286.     }
  287.   for(mptr=matrix;mptr<matrix+BLOCKSIZE;mptr++)
  288.     {
  289.       if (*mptr > 0)          /* Rounding is different for +/- coeffs */
  290.     {
  291.       *mptr = (*mptr + *qmatrix/2)/ (*qmatrix);
  292.       qmatrix++;
  293.     }
  294.       else
  295.     {
  296.       *mptr = (*mptr - *qmatrix/2)/ (*qmatrix);
  297.       qmatrix++;
  298.     }
  299.     }
  300. }
  301.  
  302. /*BFUNC
  303.  
  304. IQuantize() takes an input matrix and does an inverse quantization
  305. and puts the output int qmatrix.
  306.  
  307. EFUNC*/
  308.  
  309. void IQuantize(matrix,qmatrix)
  310.      int *matrix;
  311.      int *qmatrix;
  312. {
  313.   BEGIN("IQuantize");
  314.   int *mptr;
  315.  
  316.   if (!qmatrix)
  317.     {
  318.       WHEREAMI();
  319.       printf("No quantization matrix specified!\n");
  320.       exit(ERROR_BOUNDS);
  321.     }
  322.   for(mptr=matrix;mptr<matrix+BLOCKSIZE;mptr++)
  323.     {
  324.       *mptr = *mptr*(*qmatrix);
  325.       qmatrix++;
  326.     }
  327. }
  328.  
  329. /*BFUNC
  330.  
  331. PreshiftDctMatrix() subtracts 128 (2048) from all 64 elements of an 8x8
  332. matrix.  This results in a balanced DCT without any DC bias.
  333.  
  334. EFUNC*/
  335.  
  336.  
  337. void PreshiftDctMatrix(matrix,shift)
  338.      int *matrix;
  339.      int shift;
  340. {
  341.   BEGIN("PreshiftDctMatrix");
  342.   int *mptr;
  343.  
  344.   for(mptr=matrix;mptr<matrix+BLOCKSIZE;mptr++) {*mptr -= shift;}
  345. }
  346.  
  347. /*BFUNC
  348.  
  349. PostshiftIDctMatrix() adds 128 (2048) to all 64 elements of an 8x8 matrix.
  350. This results in strictly positive values for all pixel coefficients.
  351.  
  352. EFUNC*/
  353.  
  354. void PostshiftIDctMatrix(matrix,shift)
  355.      int *matrix;
  356.      int shift;
  357. {
  358.   BEGIN("PostshiftIDctMatrix");
  359.   int *mptr;
  360.  
  361.   for(mptr=matrix;mptr<matrix+BLOCKSIZE;mptr++) {*mptr += shift;}
  362. }
  363.  
  364. /*BFUNC
  365.  
  366. BoundDctMatrix() clips the Dct matrix such that it is no larger than
  367. a 10 (1023) bit word or 14 bit word (4095).
  368.  
  369. EFUNC*/
  370.  
  371. void BoundDctMatrix(matrix,Bound)
  372.      int *matrix;
  373.      int Bound;
  374. {
  375.   BEGIN("BoundDctMatrix");
  376.   int *mptr;
  377.  
  378.   for(mptr=matrix;mptr<matrix+BLOCKSIZE;mptr++)
  379.     {
  380.       if (*mptr+Bound < 0)
  381.     *mptr = -Bound;
  382.       else if (*mptr-Bound > 0)
  383.     *mptr = Bound;
  384.     }
  385. }
  386.  
  387. /*BFUNC
  388.  
  389. BoundIDctMatrix bounds the inverse dct matrix so that no pixel has a
  390. value greater than 255 (4095) or less than 0.
  391.  
  392. EFUNC*/
  393.  
  394. void BoundIDctMatrix(matrix,Bound)
  395.      int *matrix;
  396.      int Bound;
  397. {
  398.   BEGIN("BoundIDctMatrix");
  399.   int *mptr;
  400.  
  401.   for(mptr=matrix;mptr<matrix+BLOCKSIZE;mptr++)
  402.     {
  403.       if (*mptr < 0) {*mptr = 0;}
  404.       else if (*mptr > Bound) {*mptr = Bound;}
  405.     }
  406. }
  407.  
  408. /*BFUNC
  409.  
  410. IZigzagMatrix() performs an inverse zig-zag translation on the
  411. input imatrix and places the output in omatrix.
  412.  
  413. EFUNC*/
  414.  
  415. void IZigzagMatrix(imatrix,omatrix)
  416.      int *imatrix;
  417.      int *omatrix;
  418. {
  419.   BEGIN("IZigzagMatrix");
  420.   int *tptr;
  421.  
  422.   for(tptr=zigzag_index;tptr<zigzag_index+BLOCKSIZE;tptr++)
  423.     {
  424.       *(omatrix++) = imatrix[*tptr];
  425.     }
  426. }
  427.  
  428.  
  429. /*BFUNC
  430.  
  431. ZigzagMatrix() performs a zig-zag translation on the input imatrix
  432. and puts the output in omatrix.
  433.  
  434. EFUNC*/
  435.  
  436. void ZigzagMatrix(imatrix,omatrix)
  437.      int *imatrix;
  438.      int *omatrix;
  439. {
  440.   BEGIN("ZigzagMatrix");
  441.   int *tptr;
  442.  
  443.   for(tptr=zigzag_index;tptr<zigzag_index+BLOCKSIZE;tptr++)
  444.     {
  445.       omatrix[*tptr] = *(imatrix++);
  446.     }
  447. }
  448.  
  449.  
  450. /*BFUNC
  451.  
  452. ScaleMatrix() does a matrix scale appropriate to the old Q-factor.  It
  453. returns the matrix created.
  454.  
  455. EFUNC*/
  456.  
  457. int *ScaleMatrix(Numerator,Denominator,LongFlag,Matrix)
  458.      int Numerator;
  459.      int Denominator;
  460.      int LongFlag;
  461.      int *Matrix;
  462. {
  463.   BEGIN("ScaleMatrix");
  464.   int *Temp,*tptr;
  465.   int Limit;
  466.  
  467.   Limit=((LongFlag)?65535:255);
  468.   if (!(Temp = MakeMatrix()))
  469.     {
  470.       WHEREAMI();
  471.       printf("Cannot allocate space for new matrix.\n");
  472.       exit(ERROR_MEMORY);
  473.     }
  474.   for(tptr=Temp;tptr<Temp+BLOCKSIZE;tptr++)
  475.     {
  476.       *tptr = (*(Matrix++) * Numerator)/Denominator;
  477.       if (*tptr > Limit)
  478.     *tptr = Limit;
  479.       else if (*tptr < 1)
  480.     *tptr = 1;
  481.     }
  482.   return(Temp);
  483. }
  484.  
  485. /*BFUNC
  486.  
  487. PrintMatrix() prints an 8x8 matrix in row/column form. 
  488.  
  489. EFUNC*/
  490.  
  491. void PrintMatrix(matrix)
  492.      int *matrix;
  493. {
  494.   BEGIN("PrintMatrix");
  495.   int i,j;
  496.  
  497.   if (matrix)
  498.     {
  499.       for(i=0;i<BLOCKHEIGHT;i++)
  500.     {
  501.       for(j=0;j<BLOCKWIDTH;j++) {printf("%6d ",*(matrix++));}
  502.       printf("\n");
  503.     }
  504.     }
  505.   else {printf("Null\n");}
  506. }
  507.  
  508.  
  509. /*BFUNC
  510.  
  511. ClearMatrix() sets all the elements of a matrix to be zero.
  512.  
  513. EFUNC*/
  514.  
  515. void ClearMatrix(matrix)
  516.      int *matrix;
  517. {
  518.   BEGIN("ClearMatrix");
  519.   int *mptr;
  520.  
  521.   for(mptr=matrix;mptr<matrix+BLOCKSIZE;mptr++) {*mptr = 0;}
  522. }
  523.  
  524. /*END*/
  525.